MXM (version 0.9.4)

Orientation rules for the PC algorithm: The orientations part of the PC algorithm.

Description

The function takes the outcome of the PC algorithm, as produced by pc.skel or pc.con and performes the 4 orientation rules. A graph is also possible to visualize.

Usage

pc.or(mod, graph = FALSE)

Arguments

mod

An object with the results of the PC algorithm, as produced by pc.skel or pc.con.

graph

Boolean that indicates whether or not to generate a plot with the graph. Package RgraphViz is required.

Value

A list including: A list including:

Details

After having calculated the skeleton of the PC algorithm one may wants to perform orientations, leading to causal relationships. The rules as stated in Spirtes, Glymour and Scheines (2001) are

1. Rule 1. For each triple of vertices X, Y, Z such that the pair X, Y and the pair Y, Z are each adjacent in C but the pair X, Z are not adjacent in C, orient X - Y - Z as X -> Y <- Z if and only if Y is not in Sepset(X, Z).
2. Rule 2. If A -> B, B and C are adjacent, A and C are not adjacent, and there is no arrowhead at B, then orient B - C as B -> C.
3. Rule 3. If there is a directed path from A to B, and an edge between A and B, then orient A - B as A -> B.
4. Rule 4. This was added by Zhang (2008). If A -> B <- C, A - D - C, A and C are not adjacent, and D - B, then orient D - B as D -> B.

The first rule is applied once. Rules 2-4 are applied repeatedly until no more edges can be oriented.

References

Spirtes P., Glymour C. and Scheines R. (2001). Causation, Prediction, and Search. The MIT Press, Cambridge, MA, USA, 3nd edition.

Zhang, Jiji. (2008). On the completeness of orientation rules for causal discovery in the presence of latent confounders and selection bias. Artificial Intelligence 172(16): 1873--1896.

See Also

pc.con, pc.skel, mmhc.skel, mb

Examples

# simulate a dataset with continuous data
y = rdag(1000, 15, 0.3)
tru = y$G 
x = y$x
mod = pc.con(x)
mod$runtime

b = pc.or(mod)
plota(tru) 
dev.new()
plota(b$G)

dev.off()
plota( dag2eg(tru) )  ## essential graph
dev.new()
plota(b$G)